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ABSTRACT 

We have developed a self-consistent physical model for super stellar cluster winds 
based on combining a 1-D steady-state adiabatic wind solution and a non-equilibrium 
ionization calculation. Comparing with the case of collisional ionization equilibrium, we 
find that the non-equilibrium ionization effect is significant in the regime of a high ratio 
of energy to mass input rate and manifests in a stronger soft X-ray flux in the inner 
region of the star cluster. Implementing the model in X-ray data analysis softwares (e.g., 
XSPEC) directly facilitates comparisons with X-ray observations. Physical quantities 
such as the mass and energy input rates of stellar winds can be estimated by fitting 
observed X-ray spectra. The fitted parameters may then be compared with independent 
measurements from other wavelengths. Applying our model to the star cluster NGC 
3603, we find that the wind accounts for no more than 50% of the total "diffuse" 
emission, and the derived mass input rate and terminal velocity are comparable to 
other empirical estimates. The remaining emission most likely originate from numerous 
low-mass pre-main-sequence stellar objects. 

Subject headings: atomic processes — plasmas — stars: winds, outflows — open clusters 
and associations: individual: NGC 3603 — X-ray: general 



1. INTRODUCTION 



Super star clusters (SSCs) consist of densely-packed massive stars with typical ages of a few 
Myr. They are the most energetic coeval stellar systems, identified in a wide range of star-forming 
galaxies (Whitmore 2000). Within a characteristic radius of a few pc or less, winds from individual 
stars in such a cluster are expected to collide and merge into a so-called stellar cluster wind. In 
some extragalactic systems, such as the nuclear region of M82 (Strickland et al. 2002), winds 
from star clusters, either individually or collectively, can be energetic enough to drive galaxy-scale 
superwinds. Consequently, cluster winds can play an important role in shaping star formation and 
galaxy evolution in general (Heckman et al. 1990). 

Observationally, X-ray emission that appears to be diffuse has been detected around star 
clusters (e.g., NGC 3603, Moffat et al. 2002; Arches cluster at the Galactic center, Yusef-Zaden 
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et al. 2002, Rockefeller et al. 2005, Wang et al. 2006 ). Theoretically, both stellar cluster wind 
(e.g., Stevens & Hartwell 2003, and references therein) and galactic superwind (e.g., Strickland & 
Stevens 1999; Breitschwerdt & Schmntzlcr 1999) have been investigated, particularly in terms of 
their expected X-ray luminosities. However, neither the observed spectral nor spatial information 
has been used to confront the theoretical model prediction quantitatively. Furthermore, most 
models assume coUisional ionization equilibrium (CIE) and no non-equilibrium ionization (NEI) 
model is available for a quantitative comparison with current X-ray data of stellar cluster winds 
and/or superwinds. 

In this paper, we first introduce our non-equilibrium ionization model for stellar cluster winds 
based on combining a 1-D steady-state adiabatic wind solution (§2) with a non-equilibrium ioniza- 
tion calculation (§3). The non-equilibrium ionization effect is examined within a grid of parameter 
space, which is described in §4. Then we report its application to the diffuse X-ray emission of NGC 
3603 (§5). The physical parameters (e.g., mass input rate, terminal velocity, and metallicity) of the 
stellar cluster wind are derived from the fits and compared with other independent measurements 
and empirical estimates. 

2. X-RAY EMISSION FROM STELLAR CLUSTER WINDS 

On the scale of a galaxy, supernovae in a starburst region can drive a strong wind. Chevalier 

& Clegg (1985) describe it by an outflow model with uniformly distributed mass and energy 
depositions within a radius. Their solutions are scalable to the case of a star cluster, in which the 
mass loading and energy input are due to stellar winds instead of supernovae. 

Essentially analogous to the model of Chevalier &: Clegg (1985), Canto, Raga & Rodriguez 
(2000) present an analytic 1-D model to describe the cluster wind. In addition, their 3-D numerical 
simulation of stellar winds from 30 stars produces a mean flow that agrees well with their analytic 
model. They predict that the X-ray emission from the cluster wind of the Arches cluster (see also 
Raga et al. 2001) is detectable (Law & Yusef-Zadeh 2004). Also based on the work of Chevalier 
Clegg (1985), Stevens Sz Hartwell (2003, henceforth, SH03) present a model that allows for a 
lower energy transfer efficiency and mass loading. This model predicts global properties such as 
the X-ray luminosity and gas temperature at an individual cluster center. A crude comparison of 
the model predictions and existing Chandra observations does not give any conclusive result as to 
whether or not the diffuse emission is genuinely associated with cluster winds. 

All of the above theoretical models assume CIE for the ionization structure. However, one 
may suspect that the gas in a cluster wind may be in a highly NEI state, because of both the rapid 
shock heating in the wind-wind collision and the fast adiabatic cooling in the subsequent expansion. 
Fig. 1 illustrates such a NEI condition in a cluster wind. Clearly, the dynamic (adiabatic) cooling 
time scale defined by T/{dT/dt) is much shorter than the recombination time scales of all key ion 
species, and CIE does not hold. 
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Fig. 1. — Dynamical timescale compared with intrinsic recombination timescales for a 1-D cluster 
wind with an exponentially distributed mass input rate Mq = 2.3 x 10~^ yr~^, a terminal 
velocity Voo = 2000 km s'^, and a sonic radius = 1.38 pc. 

Moreover, the mass and energy injections in the previous theoretical models are input parame- 
ters and are crudely estimated from the properties (e.g., spectral type and radio flux) of individual 
massive cluster stars. Such a modeling approach suffers significant uncertainties and has little pre- 
dictive power. Stellar mass loss remains one of the largest uncertainties in stellar evolution (e.g., 
Kudritzki & Puis 2000). Massive stars are also not the only sources of mass loss; other cold mass 
loading, such as outflows from protostars(e.g., SH03), may also be important. 

In the following subsections, we introduce our dynamic model of a stellar cluster wind. 

2.1. Cluster Wind Dynamics 

We consider a young star cluster (age < 5 x 10^ yrs), from which the main mass and en- 
ergy injections are due to stellar winds, especially from those massive OB and Wolf-Rayet stars. 
Extending SH03's work, which assumes a uniform distribution of stars, we have constructed a 1- 
D steady-state wind model with (1) either an exponential or uniform distribution of stars; (2) a 
steady and spherically symmetric cluster wind; (3) mass and energy injections following the stellar 
mass distribution; (4) gravity and radiative cooling being neglected. This cluster wind can then be 
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described by the following steady-state Euler equations: 

1 d 



J.2 fj^j. 



{pur ) = m{r) 



du dP 

pu— = — m{r)u 

dr dr 



(1) 
(2) 



^2 (Ij. 



2^ 2^ 7 

pur \-u H — 

^ ^2 7-1/9' 



E{t) 



(3) 



where n, p and P are the velocity, density and pressure of the flow respectively, and 7 is the adiabatic 
index, taken as 5/3 for an ideal gas, while m(r) and E{r) are the mass and energy input rates per 
unit volume and are assumed to be proportional to the stellar mass density Psij^- 

The total mass and energy input rates are Jq°° A.'Kr^m.dr = Mq and Jq°° Airr^Edr = Eq. When 
the stellar mass distribution is assumed to be uniform, Mq = |7rr^m and Eq = ^irr^E, where 
Tc is the cluster radius, which coincides with the sonic radius; the analytic solutions of Eq. (1)- 
(3) can be found in Chevalier & Clcgg (1985). For an exponential stellar mass distribution, 
Ps{r) oc exp{'-r /tsc)-, where Vgc is the scale radius. We solve the above Euler equations numerically. 
Because gravity is negligible, the sonic radius depends on the stellar mass distribution alone. 

The above hydrodynamical equations do not predict how the electron temperature changes in 
the cluster. The temperature depends on how the electron is heated. We consider two extreme 
cases: one that has the electron and ion temperatures in equilibration, i.e., = Tj = T; the other 
that sets the initial T^/T = 10~^ and lets the ratio evolve via electron-ion Coulomb collisions. 
Here, the mean temperature T = [pmp/k)P/ p and p is the mean mass per particle in terms of the 
proton mass nip. Assuming cosmic abundances, the variation in the electron temperature follows 
the equation (e.g., Eq. (1) in Borkowski et al. 1994): 



d{Te/T) 
dt 



= 1.2 X 10-\{ 



T 



-3/2 



(4) 



where n is the total particle number density. 



Om model has foTir parameters: mass input rate Mq, stellar wind terminal velocity Voo (as- 
suming £^0 = 1/2MqV^), abundance Z, and the scale radius Vgc (or cluster radius for a uniformly 
mass distributed model). Fig. 2 illustrates the radial profiles of the mean/electron temperature, 
density, terminal velocity and mass input rate for these two stellar distribution models with the 
same sonic radius. Two sonic radii (0.69 pc, black; 1.38 pc, green) are considered here. When the 
sonic radius is larger, the density is lower in the inner part of the cluster, the temperature drops 
more slowly, and the terminal velocity of the wind is reached at a larger radius. 
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Fig. 2. — The radial profiles of the cluster wind (solid line — exponential distribution; dash line 
— uniform distribution) the velocity (top left panel), mean temperature (T, top right, thin lines), 
electron temperature (Tg, top right, thick lines), density (bottom left), and cumulative mass input 
rate (bottom right). Two models with different sonic radii (0.69 pc, black; 1.38 pc, green) are 
shown for each stellar mass distribution. The same mass input rate Mq = 2.3 x 10^'* Mq yr^^ and 
stellar wind terminal velocity Voo = 2000 km s^^ are assumed. 



2.2. X-ray Emission Calculation 

The broad band X-ray luminosity Lx of the cluster wind within a region (r < R) is given by: 

Lx = / A-Kr'^nenh I Kjjdv dr (5) 

where X is the frequency band of interest; Ue and are the electron and hydrogen number density 
respectively; Ky is the specific emissivity at frequency z^, and the cluster wind is assumed to be 
optically thin. Generally A^, is a function of elemental abundances, local ionic fractions and electron 
temperature. 
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The model-predicted spectrum from a given projected annulus between the inner-to-outer 
radii {Rin^Rout\ can be calculated using the following equation: 

Si,{Rin,Rout) = ^ 2TrRdR , „ ^„ — dr, (6) 

^ jRiu JR V - R^ 

where D is the distance to the cluster. The observed surface brightness can be predicted as: 

- i /; ^TPT /. L K'W.^)^^'^^^r ,7, 

where All/') and E{v' ,v) are the on-axis effective area and energy response matrix of the X-ray 
telescope / instrument . 

In principle, both the spatial and spectral information of the data can be used to constrain our 
model. But while we may reasonably assume the spectrum of the undetectable point sources to be 
the same as that of detected point sources, a similar assumption about the surface brightness cannot 
be made because the detected point sources are sparse and incomplete in a typical observation. 
Therefore, unlike the spectrum, the observed surface brightness cannot easily be compared with 
the model prediction at present. 



3. NON-EQUILIBRIUM IONIZATION SPECTRAL MODEL 

In the past, NEI models have been developed for (1) time-dependent isobaric or isochronic 
cooling of shock-heated gas (e.g., Sutherland & Dopita 1993), (2) delayed ionization in SNRs (e.g., 
Borkowski et al. 2001), and (3) delayed recombination in galactic winds (e.g., Breitschwerdt &: 
Schmutzler 1999) or in the Local Bubble (e.g., Breitschwerdt et al. 1996). The atomic data used 
are typically more than a decade old and are stored in a way that is difficult to update. Moreover, 
the physical processes included are not complete, due to either the lack of atomic data or the use 
of the CIE approximation. One example of such a process is the cascade of electrons following 
recombinations into highly excited levels, which could strongly affect the line emissivities of highly 
ionized species (e.g., Gu 2003). 

We have developed a NEI spectral model that has the following characteristics: (1) the ion- 
ization code is separated from the atomic data so that the latter can be updated conveniently; (2) 
the most updated atomic data are used; (3) recombinations into highly excited levels are included; 
(4) the electron temperature evolution, the dynamics, and the ionization structure are determined 
self-consistently. 



3.1. Atomic Data 



Our atomic data are based on CHIANTI (Young et al. 2003, version 4.2), which provides a 
database of atomic energy levels, transition wavelengths, radiative transition probabilities (A-rates) 
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and electron collisional excitation rates for 23 elements from H (Z = 1) to Zn (Z = 30), totally 206 
ions. 

Spectral calculations involving detailed fine structure levels are time-consuming and are not 
necessary for the current X-ray CCD spectral resolution. Therefore, we have combined those 
fine structure levels together according to their degeneracies and have re-calculated the collisional 
excitation rates and A-rates. For example, for the hydrogen ion there are totally 25 energy levels 
listed in CHIANTI. After regrouping, we get six principle energy levels (Is, 2s, 2p, n=3, n=4 
and n=5). For n=3, there are five sub-levels (j = 3s, 3pi/2, 3p3/2, 3^3/2, and 3d5/2), which 
we regroup as one. The new collisional de-excitation rate and A-rate are C^^i = X]j ^^'^ 
^3,« ~ X]j ^^3j,« respectively, where gj is the individual sub-level degeneracy and g is the total. 

A comparison between CHIANTI and our code with the regrouped data is shown in Fig. 3. 
The total luminosity difference between the two spectra is less than 1%. The difference in the line 
structure is very slight at > 0.5 keV and somewhat more noticeable at lower energies. But such 
difference is well within the uncertainties of the atomic data. 
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Fig. 3. — Comparison between the CIE spectra of CHIANTI (blue) and our code with regrouped 
atomic data (red), assuming T = 5 x 10^ K, cosmic abundance given by Allen (1973), and the ion 
fractions from Mazzotta et al. (1998). The black line denotes continuum. A spectrum bin size of 
lA is used. 

The ionization structure of a plasma is determined by the collisional ionization and radiative 
and dielectronic recombinations (e.g., Mazzotta et al. 1998, and references therein). CHIANTI 
does not offer the rates of these atomic processes. We have taken the collisional ionization cross- 
sections from Arnaud & Rothenflug (1985) and Arnaud & Raymond (1992), and radiative and 
dielectronic recombination rates from Verner & Ferland (1996) and Mazzotta et al. (1998), except 
for Carbon, Nitrogen and Oxygen. For these three species, we use the results of Nahar (2000), 
based on a large-scale close-coupling R-matrix calculation for photoionization cross-sections and 
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recombination coefficients. We assume tliat tfiree-body recombinations are not important at tfie 
density range under consideration and neglect them for all ions. Our current code also does not 
include charge-exchange reactions. 

Fig. 4 compares the CIE ionic fractions of our calculation with those given by Mazzotta et 
al. (1998, henceforth Mazz98). All, except for C, N, O and Fe, agree well because the same 
atomic data are used. We find good agreement for all the ions of N, with the difference being less 
than 10%. For C and O, there is good agreement except for CV, CVI, OIV, OV and OVI, for 
which the difference can reach up to 50% near the peak of maximum ionic fraction. For example, 
at T 2 X 10^ K, our ionic fraction of OV is a factor of two less than that of Mazz98. From 
T Si 3 X 10^ K to T ?s 1.5 X 10^ K where ionic fractions of Fc XV to Fe XXI are significant, our 
calculation is substantially different from Mazz98, particularly for FeXVII to FeXXII. This is due 
to different radiative recombination rates employed. 
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Fig. 4. — CIE ion fractions of our model (solid) and Mozz98's (dotted) for elements C, N, O, and 
Fe as marked in individual panels. 
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3.2. Continuum and Line Emissions 

The continuum emission consists of free-free, free-bound and two-photon continua. We include 
two important two-photon transitions: H-hke 2s Is and He-hke ls2s^So Is^ ^Sq. 

There are several ways to produce spectral lines in a plasma. Most spectral lines are produced 
by electron coUisional excitation of ions followed by the radiative decay of the excited level. At 
X-ray energies free electrons can also excite an inner-shell electron to a level above the ionization 
potential. The excited level of the ion then undergoes either a radiative decay to a lower energy 
level to produce a satellite line, or an autoionizing transition to the next ionization stage (Derc ct 
al. 2001). We follow CHIANTI in treating this autoionizing transition as a radiative decay to the 
ground level but with no emission of radiation. In addition, we include the line emission following 
an inner-shell ionization. For example, an impact ionization from the inner Is^ shell of ion 
contributes to the production of a certain line in the next higher ion We use the database of 
Mewe Sz Cronenschild (1981) for such lines (see Gorczyca et al. (2003) for an assessment of the 
database accuracy). 

Line photons are also produced during the cascade following a radiative recombination into 
an excited level. This process is particularly important in a delayed recombination scenario. A 
dielectronic recombination leads to the capture of an incident electron into an excited state above 
the ionization threshhold of the recombined ion and also emits a photon if an autoionization does 
not ensue. 



3.3. Level Balance Equations 

Electron populations at various energy levels of each ion are determined by such processes as 
collision excitation/de-excitation and radiative decay; the level balance equations are: 

rii^aij = '^njaji, (8) 

where i and j denote individual levels of an ion, aji is the total rate of the transitions between the 
two levels, and rii is the level i population relative to the total level population. 

In a low density plasma, collisional excitation generally dominates over recombination in popu- 
lating the excited states. In an over-cooled NEI state, however, the temperature may be sufficiently 
low that recombination can make a non-negligible contribution to the level population (e.g., Dere 
et al. 2001). We include radiative recombinations into the principle energy levels (n = 1 to 6) of 
each ion; the level balance equation for ion (element Z with z electrons removed) is 
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= n,(neC,, + A,,) + ne^^^^^i?i?(Z+(^+i)), (10) 

where Cji is the electron coUisional rate, Aji is the radiative decay rate (zero if j < i), and 
RR{Z^^^^^^)i is the radiative recombination rate to the level i of ion Z^^^\ Following Osterbrock 
(1989), we calculate the rate using Milne's relation and the photoionization cross sections of Karzas 
&: Latter (1961) for excited states and those of Verner &; Yakovlev (1995) for ground states. Thus, 

RR{Z+^'+% = I a,,,idu, (11) 

and 

"^"^^ - v^5.+i c2(2meA;r)3/2" ^^^^ 
where gi and (Jphot,i ^ire the degeneracy and photoionization cross section of level i of ion 



3.4. Non-equilibrium Ionization Calculations 

In general the dynamics and ionization of a plasma are entangled. The ionization structure, 
through its heating/cooling effects, can affect the temperature evolution. However, in a cluster 
wind, radiative cooling can typically be neglected and the dynamics and ionization of the wind 
can be separated. With the mass input having a certain spatial distribution in our model, and 
assuming the newly injected gas to be neutral, the cluster wind number density of the neutral ion 
Nzfi for a given element Z evolves in a cluster wind as: 

- neNz,oSz,o + neNz,iaz,i, (13) 



dr mz 

with Z and mz being the abundance and mass of element Z. 

The number density of other ions {Nz,z>o) of the same element will follow: 
1 dNz,z>o 



Up dt 



= Nz,z-iSz,z-i - Nz,z{Sz,z + az,z) + Nz,z+iOiz,z+i, (14) 



where Sz,z and az,z are the total ionization (z — > z + 1) and recombination (z — z — 1) rate 
coefficients (in cm^ s~^) of ion Our interest here is in highly ionized ions. Therefore, the exact 
ionization state of the injected gas is not important as long as it is cold or warm (i.e. T < 10^ K). 

For each element, we simultaneously solve these differential equations using IDL code LSODE 
which adaptively solves stiff and non-stiff systems of ordinary differential equations. 

In an evolving plasma, a characteristic time scale to approach a CIE state is tcie ~ ['ne{Sz,z + 
Oiz,z)]~^ (Mewe 1997). The CIE approximation requires the cooling/heating time be much longer 
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than tqie- As mentioned earlier, CIE is only a special case. If experiencing a rapid heating or 
cooling process such as thermal instability, shock, or rapid expansion, a plasma may be under- 
ionized or over-ionized, compared to the CIE result (Dopita et al. 2002). 

As an illustration, we consider a toy model for an adiabatically expanding stellar cluster wind 
with a constant mass input rate M ~ 3 x 10~^ Mq yr~^ and a constant velocity v = 1000 kms~^. 
We assume that the wind is injected at an initial radius of tq = 0.3 pc and is heated up to a 
CIE state with an equilibrium temperature Tq = 5 x 10® K. As the wind expands adiabatically, 
its temperature drops as T = To(r/ro)~^/'^. Fig. 5 compares the CIE and non-CIE spectra at 
r = 1 pc, where T ~ 1 x 10® K (left panel) and at r = 2 pc, where T ~ 4 x 10^ K (right panel). The 
difference is significant: in the lower energy part, the line intensities in the non-CIE spectrum are 
weaker than those in the CIE case, while in the upper energy part saw-tooth structures are present 
in the non-CIE spectra; all these are due to delayed recombinations. Fig. 6 shows the spectra at 




Energy (keV) Energy (keV) 

Fig. 5. — Spectra of an adiabatically expanding cluster wind at r = 1 pc (left panel) and 2 pc 
(right panel): CIE (total — blue line; continuum only — green line ) and non-CIE (total — red 
line; continuum only — yellow line) calculations. 

the same two radii with and without including the cascade of electrons following recombinations 
to n > 1 levels. The cascade lines of CV, CVI, OVII and OVIII are clearly visible. 

The more general model for a cluster wind with distributed mass and energy injections has 
been described in §2.1, and the radial profiles of the wind are shown in Fig. 2. We use the case 
of Mq = 2.3 X 10~^ Mq yr"-*^, Vqo = 2000 km s^"*^, and a sonic radius r<j = 1.38 pc to illustrate 
the effect of a distributed mass injection on the ionization structure. Fig. 7 compares the ionic 
fractions of C, N, O and Fe at four radii between the case in which Mq is all injected at the center 
and that in which it is injected in an exponential distribution. The effect of the distributed mass 
injection clearly shows up as a long ionization tail. 
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Fig. 6. — Spectra with (red) and without (blue) adding the cascades following radiative recom- 
binations to highly excited levels for an adiabatically expanding stellar cluster wind (see text for 
details). Black line denotes continuum only. 



4. MODEL PREDICTIONS FOR STELLAR CLUSTER WINDS 

To illustrate the condition under which the NEI effect is significant for stellar cluster winds, 
we here consider a cluster of an exponential stellar distribution with a scale radius r^c = 0.48 pc, 
corresponding to a sonic radius at Rs = 1.97 pc. We calculate results for a parameter grid: 
Mo = 10"^ and lO'^ M© yr'^; Kx. = 500,1000, and 2000 km s'^. The velocity, temperature, 
and density profiles of the wind are illustrated in Fig. 8. We consider both the equilibration and 
non-equilibration cases of the electron and ion temperatures (§2.1) and both the CIE and non-CIE 
scenarios for comparison. To demonstrate the results, we calculate the intrinsic cumulative spectra 
for three annular regions, [0, Rg], [Rs, ^Rg], and [2Rs, 4i?s], as well as the predicted Chandra 
ASIS-S count intensity distribution over the [0, 2Rs] range and in two energy bands, [0.3, 2.0] keV 
and [2.0, 8.0] keV. For illustration, we assume that the distance to the cluster is Z) = 10 kpc and 
that the metal abundance is solar. 

First, we find that the separate Tg evolution does not lead to significant spectral difference in 
the cumulative spectrum. It is already seen from Fig. 8 that Tg is significantly lower than T only 
when the density is low (i.e. = 2000 km s~^ and Mq = W~"^ Mq yr~^). Even in this case the 
ionic fractions of the key O and Fe ions are fairly similar (Fig. 9). The biggest difference in the 
ionic fractions between the two Tg evolutions occurs during the rise of Tg from 10^ K to 10'' K in the 
non-equilibration evolution (R < 0.01i?s). In the region making the most contribution to the X-ray 
flux {2Rs ^ R^ 0.05i?s), the ionic fractions are comparable between the two Tg evolutions in either 
the CIE or non-CIE case. Henceforth in the following discussions on the cumulative spectrum and 
count intensity distribution, only the results from the Tg evolution via Coulomb interactions are 
used. 

Second, we find that the effect of the NEI manifests itself in both the inner and outer regions 
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Fig. 7. — C, N, O and Fe ion fractions calculated for a cluster wind at four radii: panel a (0.003 pc), 
b (0.30 pc), c (0.72 pc) and d (1.38 pc). Two mass injection scenarios are assumed: all at the center 
(crosses, and the dotted lines denote the corresponding CIE calculation) and in an exponential 
distribution (triangle). Y-axis denotes the ion fraction, X-axis denotes the ionic stage (e.g., 1 
denotes I). 



of the cluster wind, but under different dynamical conditions. In the inner region it is expected 
that the CIE and non-CIE calculations will be similar for the cases of high central densities. Fig. 
10 illustrates such a case (Mo = 1 X 10"*^ Mq yr^i and Voo = 500 km s'^), in which the CIE 
approximation is good. Both calculations produce very similar cumulative spectra in either the [0, 
Rs] or [Rs, 2Rs] region, although the count intensity distributions can differ by as much as a factor 
of 2 in either the soft or hard X-ray band. On the other hand, the assumption of CIE in the inner 
region is poor for Mq = 1 x 10^'^ Mq yr^-*^ and Vqo = 2000 km s~^, which leads to a low central 
density. Fig. 11 shows clearly that the two ionization calculations produce quite different spectra. 
Also the soft X-ray count intensity from the CIE calculation is an order of magnitude weaker than 
that from the non-CIE calculation. This discrepancy is due to the difference in the ionic fraction 
of the two calculations (Fig. 9). The ionization time scale is longer than the dynamical time scale 
over which the density drops, so O and Fe are not ionized to as high a degree as specified by CIE 
at the high central temperature. We refer to this contrast as delayed ionization. Fig. 12 shows 
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Fig. 8.— Cluster wind profiles for: l^o = 500 (blue), 1000 (red), and 2000 (green) km s'^. The 
mass injection rate is the same {Mq = 10~^ Mq yr^^). The non-equilibration electron temperature 
(Te) in each case is plotted separately (dash lines). 



the spatial dependence of the ionic fraction weighted by n^R'^ pertinent to the calculation of the 
cumulative spectrum, and compares the CIE and non-CIE spectra at i? = 0.5 Rs, where the largest 
contribution to the cumulative spectrum for [0, Rs] comes from. Clearly, much stronger OVII, 
OVIII lines and Fe-L complex show up in the non-CIE spectrum. They contribute to the soft 
X-ray excess seen in Fig. 11. 

In the outer region, specifically [2Rs, 4:Rs], the density has decreased so much from spherical 
expansion that the CIE assumption is clearly not valid in any of the six wind profiles studied. The 
contrast between the non-CIE and CIE spectra then depends on how far the temperature drops 
below ~ 3 X 10^ K, when OIX ions begin to recombine. The more rapidly the temperature decreases 
with radius, the more the ionization structure deviates from the CIE result. Fig. 13 compares the 
CIE and non-CIE spectra in the [2Rs, 4:Rs] region. While the two spectra are quite similar in the 
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Voo = 2000 kms~ case, the non-CIE spectrum in the Voo = 500 kms~ case clearly produces a 
stronger hard X-ray flux from a saw-tooth-like recombination continuum. The left panel of Fig. 14 
compares Oxygen ionic fractions between CIE and non-CIE calculations. In the outer region (e.g., 
r = 3.2i?s), OIX, which produces O"'"'' recombination edges and cascade lines, is absent in the CIE 
but present in the non-CIE calculation. The right panel of Fig. 14 clearly shows the corresponding 
spectral differences in the two calculations. 

Fig. 15 shows the dependence of the cumulative X-ray luminosity on the mass input rate 
and wind velocity. The large discrepancy (up to almost one magnitude) between the CIE and 
non-CIE calculations happens for the high V^o and low Mq combination, i.e., at a high ratio of 
energy to mass input rate. In addition, the slope dependence is exactly 2 in the CIE calculation 
in which the cooling function A depends only on the local temperature and density. At a fixed 
Voo, the temperature distribution is insensitive to Mq, so L oc An^ oc Mq. In the non-CIE case, 
however, A depends on the evolutionary history of the ionization structure; the slope varies from 2 
for = 500 km s~^, in which the CIE approximation is valid, to ~ 1.3 for = 2000 km s"^. 

The above analysis indicates that the ionization structure is intimately related to the dynamical 
evolution. CIE is a fair approximation only for the case of a low ratio of energy to mass input rate. 
Some other processes not considered in our model might weaken/strengthen the non-CIE effect for 
the stellar cluster wind scenario. For example, outflows from protostars present in a cluster can 
be regarded as a mass loading process. In this case, the total mass input rate and the density will 
increase. As shown in SH03, the velocity and temperature will decrease instead, which causes a 
lower ratio of energy to mass input rate. 

Our model may also be applied to galactic superwinds. X-ray emission from strong superwinds 
may be detectable at very larger radii, where delayed recombination becomes important. But in 
this case, radiative cooling may not be neglected. Theoretical recombination spectrum has been 
shown in modelling galactic wind/outflow for the local starburst galaxy NGC3079 (Breitschwerdt 
2003). 
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Fig. 9. — Comparison of the ion fractions between the CIE (dash) and the non-CIE (sohd) calcula- 
tions of the cluster wind with Mq = 1 x 10~^ yr~^ and loo = 2000 km s~^ for = T (upper 
row) and Tg calculated via Coulomb interaction (bottom row). 
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Fig. 10. — Non-CIE model predictions (green), compared with CIE calculations (black) for the wind 
profile with Voo = 500 km s~^, as illustrated in Fig. 8. Upper panel: Cumulative spectra within 
annuli [0, Rg] (right) and [Rs, 2Rs] (left); for clarity, each cumulative spectrum from non-CIE 
calculations has been increased by a factor of 10. Lower panel: Chandra ACIS-S count intensity 
distributions in [0.3, 2.0] keV (right) and [2.0, 8.0] keV (left) bands. 
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Fig. 12. — Upper panel: Overall spectral contributions as a function of radius for Oxygen ions 
(left) and Fe ions (right) in the GIF (dash lines) and non-CIE (solid lines) scenarios. Lower panel: 
Comparison between the spectral contributions from Oxygen (left) and Fe (right) plus the total 
continuum (black) at ii = 0.5Rs- GIF (black + cyan) and non-CIF (blue+red). For clarity, the 
non-CIF spectra have been multipled by a factor of 10. All are based on the wind profile as in Fig. 
11. 
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Fig. 13. — Cumulative spectra of the cluster wind with Mq = 10^'^ Mq yr~^ and VJx) = 
500 (left panel, green) or 2000 (right panel, blue) km s~^. For comparison, the corresponding CIE 
calculations (black) are also shown. For clarity, each non-CIE spectrum has been multipled by a 
factor of 10. 
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Fig. 14. — Left panel: Comparison of the Oxygen ionic fractions between the CIE (dash) and non- 
CIE (solid) calculations of the cluster wind with Mq = 1 x 10~^ Mq yr~^ and Voo = 500 km s~^.. 
Right panel: Similar to the lower left panel of Fig. 12 but at 3.2 Rg. Free-bound recombination 
edges are denoted by the contributing ions (green). 
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Fig. 15. — CIE (dash) and non-CIE (solid) calculations of the cumulative X-ray luminosity 
versus mass input rate. An exponential stellar mass distribution with Rs = 1.97 pc and 
Voo = 500 (green), 1000 (red), and 2000 (blue) km s~^ is studied. 



- 22 - 



5. APPLICATION TO NGC 3603 

NGC 3603 at a distance of 7 kpc is the most luminous giant HII region in the Galaxy (e.g., 
Moffat et al. 1994). The star cluster responsible for the HII region is regarded as a Galactic 
analogue to SSCs observed in other galaxies. With an age of less than 3 Myr (Sung & Bessell 
2004), the cluster contains massive stars as well as low-mass pre-main sequence (PMS) stars. 

Stars in NGC 3603 have been investigated from near- infrared to X-ray (e.g.. Sung &; Bessell 
2004; Moffat et al. 2002). Based on 44 spectroscopically classified stars, which dominate the 

mass and energy injections, SH03 estimate the total mass input rate from the cluster as Mq = 
2.3 X 10^'^ Mq yr""*^ and the mean weighted stellar wind terminal velocity as Voo = 2844 km s~^. 
We use the observed diffuse X-ray emission of NGC 3603 to test our NEI model of stellar cluster 
winds. 

5.1. Data Analysis 

We use the Chandra ACIS-I observation of NGC 3603 (observation ID 633) , which was taken 
in May, 2000 for an exposure of 50 ks. Moffat et al. (2002) used this observation to study primarily 
discrete X-ray sources, although they also characterized the remaining "diffuse" emission with a 
one-temperature CIE thermal plasma model. Here, we concentrate on the diffuse emission. 

We reprocess the ACIS-I event data, using the Chandra Interactive Analysis of Observations 
software package (CIAO, version 3.2). This reprocess, including the exclusion of time intervals with 
background flares, results in a net 47.4 ks exposure (live time) for subsequent analysis. 

The superb spatial resolution of the data allows for the detection and removal of much of the X- 
ray contribution from individual stars. Wc search for X-ray sources in the 0.5-2, 2-8, and 0.5-8 kcV 
bands by using a combination of three source detection and analysis algorithms: wavelet, sliding- 
box, and maximum likelihood centroid fitting (Wang 2004). Source count rates are estimated with 
the 90% energy-encircled radius (EER) of the point spread function (PSF). The accepted sources 
all have the local false detection probability P < 10~^. The right panel of Fig. 16 shows the 
detected point sources within the cluster field. The source detection in the inner ~ 20" radius of 
the cluster is highly incomplete due to severe source confusion. 

To analyze the remaining "diffuse" emission from the cluster, we remove the detected source 
regions from the data. For each faint source with a count rate (CR) < 0.01 counts s~^, we exclude 
a region of twice the 90% EER, which increases with off-axis angle. For each source with CR 
> 0.01 counts s~^, an additional factor of 1 -|- log(CR/0.01) is multiplied to the source-removal 
radius to further minimize the residual contribution from the PSF wings of relatively bright sources. 

We extract the spectra of the diffuse X-ray emission separately from the two annular regions, 
Da and Db (divided to contain roughly the same number of counts), as well as an accumulated point 



-23- 




RIGHT ASCENSION (J2000) 

Fig. 16.— Chandra ACIS-I observation of NGC 3603. Left panel: The 0.5-8.0 keV intensity 
contours (at 52, 58, 70, 94, 130, 178, 238, 310, 394, and 490 xlO""^ counts s"^ arcmin-^) overlaid 
on an optical image. Right panel: A close-up of the central cluster. The image shows the ACIS-I 
event distribution. Point-like X-ray sources are marked individually with small circles of 2x the 
PSF 90% EER. The two annulii. Da ([19", 30"]) and Db ([30", 83"]), are used to extract the 
spectra of the diffuse emission. 

source spectrum from the combined region (Da+Db). An additional spectrum extracted from the 
outer annulus between A'.O and 5'0 radii is used as the local background of the diffuse X-ray emission, 
which in turn are normalized for the background subtraction of the sources. Fig. 17 shows a joint- 
fit of the spectra with our cluster wind model plus a thermal model, which is used to characterize 
the point-like source spectrum. Here wc assume that the detected and undetected point sources 
have the same spectral shape. We jointly fit the diffuse spectra and the accumulative spectrum 
of point sources together. Wc find that the spectral fit results arc insensitive to the assumed 
stellar distributions of our model and to the exact scale radius. Therefore, Table 1 includes only 
the results from the model with the more physical exponential stellar mass distribution; the fixed 
scale radius Vgc = 0.68 pc corresponds to the same sonic radius as adopted by SH03. The best-fit 
model predicts the luminosity contributions of L2-S keV ~ 2.5 x 10^^ erg s~^ and 8.8 x 10^^ erg s^^ 
from the cluster wind and residual point-like source components. Fig. 18 shows the predicted 
radial surface intensity distribution of the cluster wind contributions compared with the observed 
"diffuse" emission in two energy bands. 
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Table 1. Joint spectral fit results^ 



Parameter Model 





(10-- cm--) 


81+"-"^' 




Cluster wind component 




(pc) 


0.68 (fixed) 


Abundance 


"-'•^-'■-0.11 


LogMo 


(Mo yr) 


07+0.42 




(1 X 10^ km s"^) 


1.53 (> 0.50) 




Point-like source component 


kT 


(keV) 


7Q+0.41 

^•'^-0.33 




(10~^ photons cm^^) 


«Q+0.77 
''•'-'"J-1.56 


K2' 


(10~^ photons cm^^) 


' - '0-1.73 




(10~^ photons cm~^) 


9/^+0.22 
"^•^"-0.19 


xVd.o.f. 


301.3/457 



^Note: all error bars are at the 90% confidence 

^Ki is the MEKAL normalization of the undetected 
point sources in Region Da. 

is same as Ki, but in Region Db. 

'^K^ refers to the MEKAL normalization of the de- 
tected point sources. 
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Fig. 17. — Joint fit to the spectra of the point-like sources (green) and the two diffuse X-ray spectra 
(red and blue) of NGC 3603. All these spectra are grouped to contain at least 25 net counts per bin. 
The point-source spectrum is modelled with an optically-thin thermal plasma, while the diffuse 
X-ray spectral model also contains a joint-fitted residual point-like source contribution (dotted 
histograms) as well as the cluster wind component (dashed). 

5.2. Discussion 

The above results represent the first quantitative comparison of a self-consistent NET cluster 
wind model with observational data. Our fitted foreground X-ray-absorbing gas column density, 
Nh = 8.lt°oi X 1021 cm 2, is consistent with the interstellar reddening toward NGC 3603. The 
reddening of E{B — V) = 1.25 — 1.8 mag, determined by Sung & Bessell (2004) from 115 bright 
stars in the field, corresponds to Nh = (7.3 - 10.4) x 10^^ cm"^ (e.g., Bohlin et al. 1978). The 
derived abundance of the wind is sub-solar: Z = 0.21^q'}^ Zq. This abundance likely represents 
an underestimate of the true value due to the oversimplification in the 1-D model, which does not 
account for substructures in both temperature and density. Such substructures lead to smearing 
out of spectral line features that characterize the metal abundance. Our fitted mass input rate 
and terminal velocity of the cluster wind is, within a factor of ~ 2, consistent with the empirical 
estimation by SH03 (§5). But the cluster wind contributes only about 25% of the observed total 
diffuse emission luminosity. 
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Fig. 18. — ACIS-I intensity profiles of the diffuse X-ray emission (crosses with la error bars) in 
the 1.0-2.5 keV (left) and 2.5-4.0 keV (right) bands. The solid lines represent the predicted cluster 
wind contributions based on the spectral fit. 

While the remaining of the diffuse emission likely originate from undetected point-like sources, 
what is their nature? Moffat et al. (2002) and Sung & Bessell (2004) suggested that PMS young 
stellar objects (YSOs) may be responsible. According to the study by Saga et al. (2001), the 
NGC 3603 cluster has an initial mass function (IMF) that can be characterized by a power law 
with a slope of 0.84 and contains about 2 x 10^ stars in the mass range 7-75 Mq. Using this IMF, 
we estimate the presence of about 8 x 10^ YSOs (0.3-3 M0) in this cluster. The study of such 
YSOs in the Orion Nebular (Feigelson et al. 2005) shows that their average 2-8 kcV luminosity is 
~ 1 X lO^'^erg s~^ per star. Thus the expected total YSO contribution in the NGC 3603 cluster is 
L2-8 keV '"^ 8 X 10^^, consistent with the prediction from our spectral fit. 

6. SUMMARY 

In this paper, we have presented a NEI spectral code, which is based on our updated atomic 
data, including the atomic process of recombinations into highly excited levels. We have shown 
important differences between the CIE and non-CIE calculations. 

We have constructed a self-consistent physical model for super stellar cluster winds, combining 
the NEI code with a 1-D steady-state adiabatic wind solution. We find that in the inner region of 
the star cluster, the NEI effect is significant in the regime of a high energy to mass input ratio and 
manifests in an enhanced soft X-ray emission, which could be an order of magnitude more luminous 
than that obtained under the CIE assumption. 

We apply the cluster wind model to the Chandra data on NGC 3603 and find that the inferred 

mass and energy input rates of the stellar cluster winds are comparable to the empirical estimates. 
But the cluster wind in NGC 3603 accounts for no larger than 50% of the total diffuse X-ray 
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emission. YSOs are likely responsible for the remaining diffuse emission. 

We are grateful to Z. Li for providing the software programs used for calculating the wind 
dynamics. We thank K. J. Bowkowski, M.-F. Gu, J. S. Kaastra, F. Paerels, and R. Smith for many 
valuable comments and suggestions on this work, which is supported by NASA through the grant 
SAO/CXC GO4-5010X. 
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